Molecular-spin dynamics study of electromagnons in multiferroic RMn205 
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We investigate the electromagnon in magnetoferroelectrics RMn2 05 using combined molecular- 
spin dynamics simulations. We identify the origin of the electromagnon modes observed in the 
optical spectra and reproduce the most salient features of the electromagnon in these compounds. 
We find that the electromagnon frequencies are very sensitive to the magnetic wave vector along 
the a direction. We further investigate the electromagnon in magnetic field. Although the modes 
frequencies change significant under magnetic field, the static dielectric constant electromagnon does 
not change much in the magnetic field. 
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Magnetoferrelectrics with strong magnetoelectric 
(ME) couplings have attracted intensive attention for 
their novel physics [l|-y| and potential applications in 
tunable multifunctional devices [5|, |6|. One of the most 
attractive features of the magnetoferroelectrics is that 
the spin waves and infrared phonons could couple via 
ME interactions leading to a new type of elementary 
excitation, electromagnon [7[ . Electromagnon has mixed 
characters of both magnons and infrared phonons. It 
therefore opens a way to control the magnetic excita- 
tions via electric fields [8*, which may have important 
applications for information processes. However, despite 
electromagnon has been proposed theoretically almost 
four decades ago, it has been observed experimentally 
only very recently [9J. Till today, our knowledge about 
electromagnon is still very limited. 

In the past few years, electromagnon has been observed 
in various magnetoferroelectrics, mostly in the family of 
RMnOaji, [l^ and RMnaOg [l|] (R =Tb, Dy, Gd, Eu, 
Y) compounds. So far, most theoretical works on electro- 
magnon were concentrated on the RMnOs family |12l4l4| . 
because they have relatively simpler magnetic structures, 
containing only Mn'^^ ions. The study of electro mag non 
in RMn205 is still very few both experimentally ll|, |l5| 



and theoretically [16[ , due to the complexity of the mag- 
netic structures in these compounds. A simulation of the 
electromagnon in RMn205 at atomistic level is still lack. 
It has been shown experimentally that the dielectric step 
at the commensurate (CM) to incommensurate (ICM) 
magnetic phase transition in RMn205 is caused by the 
electromagnon ^11.] . However, it is still unclear how the 
electromagnon is related to the magnetodielectric effect, 
which also happens at the same temperature and is one 
of the most prominent effects in these materials [Ij . 

In this work, we carry out combined molecular-spin 
dynamics (MSD) simulations to study electromagnons in 
RMn205. We reproduce the most salient features of the 
electromagnon in these compounds and find that the elec- 
tromagnons frequencies are very sensitive to the magnetic 
wave vector along the a direction. We further investigate 
the behavior of electromagnon under magnetic fields, and 
find that although the mode frequency change significant 
under magnetic field, the static dielectric constant elec- 
tromagnon does not change much in the magnetic field 
in our model. 



To study the electromagnon in RMn^O 
effective Hamiltonian derived in Ref. 
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where J3 is the Mn''+- Mn^+ superexchange interac- 
tion through pyramidal base corners, and J4 the su- 
perexchange interaction through the pyramidal apex [19j . 
Mn'^+ couples to Mn^+ either antiferromagnetically via 
J4 along a axis or with alternating sign via J3 along 
h axis, whereas Mn^+ ions in two connected pyramids 
couple each other antiferromagnetically through J5. Ji, 



J2 couple Mn^+ ions along the c axis, u is an artifi- 
cial optical phonon to describe the electric polarization 



171 |18|. In RMn2 05, the infrared modes that couple 



to the electromagnon are the same modes that causes 
the marcoscopic polarization. This is very different from 
those in RMnOs compounds. One can find more details 
of the model in Ref. Ig]. The major difference between 



the current raodel and previous ones[17|, [iSJ is that we 
include here the single ion easy axis anisotropy K^ for 
Mn'^+ ions and the easy plane anisotropies Di for Mn'^+ 
and Mn^"*" ions to describe the magnetic anisotropy in 
these compounds. H is the applied magnetic field and g 
is gyromagnetic ratio. 

We use the model parameters obtained in Ref. [l8| un- 
til otherwise noticed. Since accurate magnetic anisotropy 
energies are still very hard to obtain from the current 
first-principles methods, especially for complex materials 
like RMn205, we use empirical anisotropy parameters to 
fit the experimental magnetization and optical spectra. 

We carry out MSD simulations of the Hamiltonian 
by solving the coupled equations of motion of spin and 
phonon using forth-order Runge-Kutta method, 



m, - -s, X Hf ^ 

ilk = Ffc/m 



(2) 
(3) 



where 'H-l — —dH/dSi is the effective magnetic field 
acting on the zth Mn spin S^ and Ffc = ~dH/duk is 
the force acting on the kth local phonon mode. Since 
Mn'^+ and Mn'*+ have 4 and 3 unpaired local d elec- 
tron respectively, we set the norm of MnjJ" spin vector 
I S3 1 — -\/2(2 -I- 1) and the norm of MnJ spin vector 
IS4I = 1/1.5(1.5 -|- 1). Unlike previous methods JJ,], our 
method takes the phonon degree of freedom into consid- 
eration explicitly. 

We start the MSD simulation from a pool of equilib- 
rium configurations sampled from MC simulation at a 
given temperature T. The system is further relaxed for 
sufficient time in the MSD simulations. Simulations with 
different initial configurations are averaged to realize a 
canonical assemble at fixed T [20] . The optical spectrum 
Imx(w) is calculated from the Fourier-Laplace transform 
of correlation function $(i) of dipole moment P(t) with 
<I>(t)^(P(t)P(0))-(P(i))(P(0)),i.e., 



Imxi'^) 



KbT 



(4) 



where P{t) = ^j, Uk{t) and (• • • ) represents the ensemble 
average. The static susceptibility x(0) can be calculated 
as, 



X(0) 



KbT' 



(5) 



where Kb is the Boltzman constant. To extract the de- 
tailed parameters of the spectral weight, we fit x(w) to a 
lorentz model, 



j J 'J 



(6) 



where j enumerates the vibrations, Wj is the resonance 
frequency, and 7^ is the damping rate. To avoid numer- 
ical errors, we fit directly the correlation function $(i). 
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FIG. 1: (Color online). The imagery part of the dielectric 
response Imx(tj) calculated using the 12 x 12 x 12 lattice in 
the CM magnetic phase (red) and the ICM magnetic phase 
(black). Inset: The electromagnon frequency as function 
of the magnetic wave vector along a direction q^. Circles 
represents the electromagnon frequency calculated by MSD, 
whereas the black line is calculated from J5 sin(27rga:) — 14.63 



which can be analytically Fourier transformed to lorentz 
model in Eq. [51 The static electric susceptibility is the 
sum of all oscillators strengths, x(0) = X)? ^- Accord- 
ing to the sum rule, the total spectral weight 5* = ^ Sj 
should be conserved. 

We first compare the MSD results to previous MC sim- 
ulations [18|, by carrying out the MSD simulation in a 
12x12x12 lattice at various temperatures, without in- 
cluding the magnetic anisotropy. The imaginary part of 
the dielectric response Imx(a;) of the CM phase at 20K 
and ICM phase at 3K is shown in Fig. [T] We can see 
that, in the CM phase, there is only one peak centered 
at about 185 cm~^ , corresponding to the bare phonon fre- 
quency in our model. While at ICM phase, there are two 
peaks: one is the phonon peak at about 193 cm~^ which 
is hardened compared to the bare phonon frequency due 
to strong spin-phonon coupling [21j, and the new peak 
appear at about 58 cm~^ is the electromagnon. 

The calculated x(0) is shown in Fig. HJa) at different 
temperatures compared with those obtained from MC 
simulations |18|. We see that x(0) calculated from both 
MC and MSD shows a step at the CM to ICM transi- 
tion. The difference between x(0) calculated with MC 
an MSD comes from the long-time relaxation effects that 
are not included in the MSD simulation. The results are 
in remarkable agreement with the experimental results 
of Ref. [llj. Although electromagnons have been suc- 
cessfully reproduced using our model, the frequency of 
the electromagnon mode is much higher than the exper- 
imental values for TbMn205 and YMn9 0f^[lH. Further 
examinations show that the electromagnon frequencies 
are very sensitive to the x component of wave vector qx . 
In the ICM phase, the spins connected by the strongest 
interaction J4 are almost parallel to each other [18j, and 
the spiral along a direction comes from the spins con- 
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FIG. 2: (Color online) (a) The static dielectric susceptibility 
x(0) as functions of temperature calculated using MC and 
MSD. (a) The imaginary part of the dielectric response Im 
e{ij) of RMn205 under different magnetic fields. 



nected by J5, with an angle of iirq^. Since the electro- 
inagnon excitation is at Brillouin zone center, J4 makes 
no contribution to the electromagnon as can be easily 
seen from Eq. [2l whereas the second largest exchange 
interaction J5 dominates the interaction. We therefore 
expect the energy of the electromagnon ~ Jc,sa\{2'Kqx). 
To confirm this result, we simulate the electromagnon 
on a larger 50x6x6 lattice and adjust the parameter 
slightly to obtain different qx of the ICM state. The 
calculated lowest electromagnon mode frequencies at dif- 
ferent qx are shown in the insert of Fig. [1] As we 
see that the electromagnon frequencies decrease signif- 
icantly as qx approaching 7r/2, showing nice agreement 
with J5sin(27rq^) — 14.63cm^^. At qx=OA%, the electro- 
magnon frequency decreases to about 10 cm~^, which is 
in good agreement with experiments [11[. Most of the 
RMn205 (R=Tb, Ho, ...) compounds have qx ~ 0.48, 
therefore the electromagnon frequencies in these mate- 
rials should be very similar. However, we expect that 
the electromagnons in TmMn2 05 24], with q^ ~ 0.467, 
should have much higher frequencies. Experimental con- 
firm of this prediction is called for. 

We now look into more details of the electromagonons 
in the ICM phase. To compare better with the ex- 
periments, we slightly modify our model parameters to 
J3=1.5 meV and J3=0.3 meV to get ga;=0.46 in the 
50x6x6 lattice [23[. Though the frequencies are a little 
higher than the experimental data, it does not change 
the basic physics of our model. To see how the mag- 
netic anisotropy affect the electromagnon, we also include 
the ionic anisotropy energies in the simulation. We take 
D3 = 0.09 meV, Da = 0.06 meV and three typical in- 
plane easy axis anisotropy energies X3=(0.19, ± 0.21, 
0), (0.278, ±0.082, 0), (0.08 ±0.28, 0), corresponding to 
the case of Dy, Tb and Ho respectively. We first show 
the results using ii:3=(0.19, ± 0.21, 0). 

The imaginary part of the dielectric response Im e{uj) 
at 3K are presented in Fig. [21 After adding the mag- 




FIG. 3: (Color online) The cartoons of two identified electro- 
magnon modes. The blue arrows denote rotation direction. 
(a) The "phason" mode, where spins in the same chain rotate 
in the same direction as a unit, (b) The highest frequency 
electromagnon mode appear under magnetic field H || a. The 
blue arrows denote rotation directions. 



netic anisotropy, the original electromagnon peak splits 
into three modes , centered at 30 cm~^, 46 cm~^, 51 
cm~\ respectively, consistent with recent experiments 
[nl [l5| . To extract the electromagnon modes, a short 
electrical pulse E \\ & at i = is applied at the equilib- 
rium ground states. The time evolution of each spin is 
calculated and stored. The electromagnon normal modes 
are reconstructed from the Fourier transform of each spin 
motion. From the calculated mode profiles, we confirm 
that the strongest and lowest frequency peak is the opti- 
cal counterpart of "phason" with wave vector Q=0. The 
pattern of the mode is shown in Figl^fa), which corre- 
sponds to the spin rotation in the spiral plane [15| . In 
this mode, the spins of neighboring AFM chains rotate in 
opposite directions, and couple to phonon modes through 
the term J3 mS^ • Sj forming an electromagnon. We do 
not obtain a clear profile of the other two modes, because 
they are relatively weak and mixed with other modes. 

Magnetodielectric effect is one of the most prominent 
effects in RMn205. To clarify the relation between the 
electromagnon and magnetodielectric effect, it is impor- 
tant to study the electromagnon in the magnetic field. 
Figured] depicts the imaginary part of dielectric response 
for H=8, 12, 16 T II a direction. Interestingly, there 
is a new electromagnon peak emerges in the magnetic 
field above the three electromagnon modes of the zero 
magnetic field. The new mode can be characterized as 
spins rotating around an axis perpendicular to itself in 
the ab plane, as shown in Fig|3Ub). We find that the pat- 



tern of the lowest frequency mode does not change, and 
the magnitude of the intermediate two modes are getting 
stronger in magnetic field. The frequencies of the two 
lowest electromagnon peaks shift down with the increase 
of magnetic field and the frequency of the lowest mode 
shifts down by about 17% under H^16T compared to 
that of zero field. In contrast, the two high frequency 
modes shift toward even higher frequency with increas- 
ing of the magnetic field. To identify the driving force 
for electromagnon frequency shifts, we compare the mag- 
netic structure with and without magnetic field and find 
that they are almost the same. Therefore, the frequency 
shifts mainly come from the dynamic contribution of Zee- 
man term —g ^^ H • S^. 

In spite of the large frequency change in magnetic field, 
the total static susceptibility calculated from MSD re- 
main almost unchanged compared with that at H=0. 
The total spectral weight also conserves well as the mag- 
netic field increases, within numerical errors. The most 
contribution to x(0) always comes from the lowest elec- 
tromagnon. At zero magnetic field, the intermediate two 
electromagnon modes have negligible contribution to the 
overall x(0)- Under magnetic field, the lowest electro- 
magnon transfers spectral weight to higher ones, e.g. 50% 
at H=16 T. Therefore although the frequency decreases 
significantly, x(0) calculated from MSD almost keep con- 
stant. 

We also study the electromagnon in magnetic fields 
with H|| 5 and H || c. The results for H || 5 are similar to 
those for H || a, while H j| c has weak influence on elec- 
tromagnon both in frequency and spectral weight. Our 
simulations suggest that the magnetodielectric effects in 
TbMn205 are not caused by the electromagnon, but it 
does not rule out such possibility for DyMn205, which 
has unique magnetic structure in the ICM phase, with 
qx—0.5. We also do the simulations using two other K3 
parameters, and obtained very similar results. 

The behaviors of electromagnons in RMn2 05 are very 
different from those in DyMnOs where the electromagnon 
shows a soft-mode behavior in magnetic field H \\ a 
and lead to the increase of dielectric constant [24| in the 
field. This is probably because the magnetic structures 
of DyMnOa are much more sensitive to the applied mag- 
netic field. 

To summarize, we investigated the electromagnon 
in magnetoferroelectrics RMn205 using combined 
molecular-spin dynamics simulations. We have identi- 
fied the origin of the electromagnons in these compounds, 
and reproduced the most salient features of the electro- 
magnon in these compounds. We find that the electro- 
magnon frequencies are very sensitive to the magnetic 
wave vector along the a direction. We further investi- 
gate the electromagnon in magnetic field. Although the 
mode frequency change significant under magnetic field, 
the static dielectric constant does not change much in the 
magnetic field in our model. 
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